Autocorrelation spatiale

Application aux résultats des élections européennes de 2024 en France

Auteur·rice

Claude Grasland

Date de publication

2024-12-11

Vous allez choisir une des listes électorales puis reprendre l’ensemble des analyses vues en cours à l’aide des programmes ci-dessous. Les choix possibles sont :

“Aubry” “Bardella” “Bellamy” “Deffontaines” “Glucksman” “Hayer” “Lassale” “Marechal” “Toussaint” “Autres”

varY <- "Bardella"
nameY<-paste("% vote", varY)

A. ECHELLE REGIONALE

Préparation des données

### Fonds de carte
map<-st_read("data/elect2024/map_reg.shp", quiet=T) %>%
      select(code=code_reg,  geometry) %>%
      st_transform(2154)

### Données
don<-readRDS("data/elect2024/don_reg.RDS") %>%
          select(code=code_reg, nom=nom_reg, 3:12)
don$nom<-as.factor(don$nom) 
levels(don$nom)<-c("ACAL", "AQUI","AURA","BOFC","BRET","CVDL","IDF","OCCI","NOPI","NORM","PDL","PACA")

### Jointure
mapdon<-left_join(map[,1],don)

# Table de voisinage
map_nb<-spdep::poly2nb(mapdon,row.names = mapdon$code)

# Table de poids
map_nb_w<-nb2listw(map_nb)
summary(map_nb_w)
Characteristics of weights list object:
Neighbour list object:
Number of regions: 12 
Number of nonzero links: 46 
Percentage nonzero weights: 31.94444 
Average number of links: 3.833333 
Link number distribution:

2 3 4 5 6 
2 3 3 3 1 
2 least connected regions:
53 93 with 2 links
1 most connected region:
24 with 6 links

Weights style: W 
Weights constants summary:
   n  nn S0       S1       S2
W 12 144 12 6.643889 49.32222
# Carte de voisinage
coo<-st_coordinates(st_centroid(mapdon))
mf_map(mapdon, type="base",col="lightyellow")
mf_layout("Carte des liens de voisinage", frame=T)
plot.nb(map_nb,coords = coo,add = T,col = "red",points = F)
mf_label(mapdon, var="code", col="blue",cex=0.6,halo = T, bg="white",r = 0.1)

Choix d’une variable Y

sel<-mapdon %>% select(code,nom, Y=varY) %>%
                mutate(Y_std = as.numeric(scale(Y)) ) 
par(mfrow=c(1,2))
# Carte de Y
mf_map(sel, var="Y", type="choro",
       nbreaks = 6,
       leg_title = "Quantiles",
       leg_val_rnd = 1)
mf_layout("Y",
          frame=T,
          credits="",
          arrow=F)

# Carte de Ystd
mypal<- brewer.pal(6, "RdYlBu")
mybreaks<-c(-10,-2,-1,0,1,2,10)
mf_map(sel, var="Y_std", type="choro",
       breaks=mybreaks, 
       pal = mypal,
       leg_title = "moy. & ec.type",
       leg_val_rnd = 1
       )
mf_layout("Y_std",
          frame=T,
          credits="",
          arrow=F)

Valeurs de voisinage

sel$Y_lag<-as.numeric(lag.listw(map_nb_w,sel$Y))
sel$Y_std_lag<-as.numeric(lag.listw(map_nb_w,sel$Y_std))
#mapsel$Y_std_lag<-lag.listw(map_nb_w,mapsel$Y_std)
st_drop_geometry(sel)
   code  nom        Y      Y_std    Y_lag   Y_std_lag
1    11  IDF 18.78888 -2.1535547 37.61953  0.72909500
2    24 CVDL 34.93906  0.3187601 30.11507 -0.41970883
3    27 BOFC 37.08970  0.6479876 30.74364 -0.32348574
4    28 NORM 35.32945  0.3785223 29.86869 -0.45742619
5    32 NOPI 42.40672  1.4619321 30.81702 -0.31225318
6    44 ACAL 38.33273  0.8382729 32.76177 -0.01454500
7    52  PDL 27.62945 -0.8002157 31.69672 -0.17758538
8    53 BRET 25.57934 -1.1140528 31.47945 -0.21084671
9    75 AQUI 30.93906 -0.2935711 31.79258 -0.16291148
10   76 OCCI 33.68791  0.1272309 33.49939  0.09837241
11   84 AURA 30.91391 -0.2974213 35.06019  0.33730344
12   93 PACA 38.64522  0.8861097 32.30091 -0.08509518

Indice et diagramme de Moran

spdep::moran.test(sel$Y,map_nb_w, alternative = "two.sided")

    Moran I test under randomisation

data:  sel$Y  
weights: map_nb_w    

Moran I statistic standard deviate = -0.69097, p-value = 0.4896
alternative hypothesis: two.sided
sample estimates:
Moran I statistic       Expectation          Variance 
      -0.20852369       -0.09090909        0.02897351 
par(mfrow=c(1,1))
moran.plot(x=sel$Y_std,
           listw=map_nb_w, 
           pch=20,
           main ="Diagramme de Moran",
           xlab = "Variable observée (Y)", 
           ylab=  "Variable de voisinage (Ylag)",labels = F)
text(sel$Y_std, 
     sel$Y_std_lag, 
     sel$nom, 
     cex=0.5, 
     col="red",
     pos=1)

Calcul des LISA

locm<-localmoran(sel$Y,map_nb_w,alternative = "two.sided")
sel2<-as.data.frame(locm)
tabres<-cbind(sel,sel2)
st_drop_geometry(tabres)
   code  nom        Y      Y_std    Y_lag   Y_std_lag          Ii         E.Ii
1    11  IDF 18.78888 -2.1535547 37.61953  0.72909500 -1.71288650 -0.459946894
2    24 CVDL 34.93906  0.3187601 30.11507 -0.41970883 -0.14594885 -0.010076829
3    27 BOFC 37.08970  0.6479876 30.74364 -0.32348574 -0.22867063 -0.041641779
4    28 NORM 35.32945  0.3785223 29.86869 -0.45742619 -0.18888653 -0.014209497
5    32 NOPI 42.40672  1.4619321 30.81702 -0.31225318 -0.49799230 -0.211958224
6    44 ACAL 38.33273  0.8382729 32.76177 -0.01454500 -0.01330111 -0.069689401
7    52  PDL 27.62945 -0.8002157 31.69672 -0.17758538  0.15502539 -0.063505302
8    53 BRET 25.57934 -1.1140528 31.47945 -0.21084671  0.25624840 -0.123085644
9    75 AQUI 30.93906 -0.2935711 31.79258 -0.16291148  0.05217394 -0.008547175
10   76 OCCI 33.68791  0.1272309 33.49939  0.09837241  0.01365383 -0.001605393
11   84 AURA 30.91391 -0.2974213 35.06019  0.33730344 -0.10944134 -0.008772835
12   93 PACA 38.64522  0.8861097 32.30091 -0.08509518 -0.08225854 -0.077870117
        Var.Ii        Z.Ii Pr.z....E.Ii..
1  0.357689878 -2.09496520     0.03617409
2  0.009975286 -1.36040224     0.17370267
3  0.083806257 -0.64605635     0.51824287
4  0.020170926 -1.22990878     0.21873126
5  0.534502193 -0.39123980     0.69561999
6  0.207464924  0.12379885     0.90147454
7  0.124891996  0.61836534     0.53633454
8  0.582852069  0.49686985     0.61928083
9  0.017795654  0.45517954     0.64898008
10 0.005129009  0.21306678     0.83127487
11 0.012522057 -0.89961308     0.36832619
12 0.387754354 -0.00704742     0.99437702

Carte des LISA

par(mfrow=c(1,2))

# Carte de Ystd
mypal<- brewer.pal(6, "RdYlBu")
mybreaks<-c(-10,-2,-1,0,1,2,10)
mf_map(sel, var="Y_std", type="choro",
       breaks=mybreaks, 
       pal = mypal,
       leg_title = "moy. & ec.type",
       leg_val_rnd = 1
       )
mf_layout("Variable standardisée",
          frame=T,
          credits="",
          arrow=F)


mypal<-brewer.pal(5,"BrBG")
mybreaks<-c(-10,-2,-1,1,2,10)
mf_map(tabres, var="Z.Ii", type="choro",
       breaks=mybreaks, 
       pal = mypal,
       leg_title = "Z-score",
       leg_val_rnd = 1
       )
mf_layout("LISA",
          frame=T,
          credits="",
          arrow=F)

Typologie de Moran

# create a new variable identifying the moran plot quadrant for each observation, dismissing the non-significant ones
tabres$typ <- NA

# high-high quadrant
tabres[(tabres$Y_std >= 0 & 
                 tabres$Y_std_lag >= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "high-high"
# low-low quadrant
tabres[(tabres$Y_std <= 0 & 
                 tabres$Y_std_lag <= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "low-low"
# high-low quadrant
tabres[(tabres$Y_std >= 0 & 
                 tabres$Y_std_lag <= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "high-low"
# low-high quadrant
tabres[(tabres$Y_std <= 0 & 
                 tabres$Y_std_lag >= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "low-high"
# non-significant observations
tabres[tabres$Pr.z....E.Ii.. > 0.05, "typ"] <- "not signif."  

tabres$typ <- as.factor(tabres$typ)
tabres$typcol<-NA
tabres$typcol <- case_when(tabres$typ=="high-high" ~ "red",
                        tabres$typ=="low-low" ~ "blue",
                        tabres$typ=="high-low" ~ "pink",
                        tabres$typ=="low-high" ~ "lightblue",
                        tabres$typ=="not signif." ~ "gray90"                      
                        )

mf_map(tabres, type="base",  col=tabres$typcol)
mf_layout("Typologie de Moran", 
          frame=T,
          credits ="")
mf_legend(type="typo",val = c("haut-haut","bas-bas","haut-bas","bas-haut","non sign."),pal = c("red","blue","pink","lightblue","gray90"),title = "Types")

B. ECHELLE DEPARTEMENTALE

Préparation des données

### Fonds de carte
map<-st_read("data/elect2024/map_dept.shp", quiet=T) %>%
      select(code=code_dpt,  geometry) %>%
      st_transform(2154)

### Données
don<-readRDS("data/elect2024/don_dept.RDS") %>%
          select(code=code_dpt, nom=nom_dpt, 5:14)


### Jointure
mapdon<-left_join(map[,1],don)

# Table de voisinage
map_nb<-spdep::poly2nb(mapdon,row.names = mapdon$code)

# Table de poids
map_nb_w<-nb2listw(map_nb)
summary(map_nb_w)
Characteristics of weights list object:
Neighbour list object:
Number of regions: 94 
Number of nonzero links: 474 
Percentage nonzero weights: 5.364418 
Average number of links: 5.042553 
Link number distribution:

 2  3  4  5  6  7  8 10 
 6 12 14 19 30 11  1  1 
6 least connected regions:
06 29 57 62 66 74 with 2 links
1 most connected region:
77 with 10 links

Weights style: W 
Weights constants summary:
   n   nn S0       S1       S2
W 94 8836 94 40.29206 383.9874
# Carte de voisinage
coo<-st_coordinates(st_centroid(mapdon))
mf_map(mapdon, type="base",col="lightyellow")
mf_layout("Carte des liens de voisinage", frame=T)
plot.nb(map_nb,coords = coo,add = T,col = "red",points = F)
mf_label(mapdon, var="code", col="blue",cex=0.6,halo = T, bg="white",r = 0.1)

Choix d’une variable Y

sel<-mapdon %>% select(code,nom, Y=varY) %>%
                mutate(Y_std = as.numeric(scale(Y)) ) 
par(mfrow=c(1,2))
# Carte de Y
mf_map(sel, var="Y", type="choro",
       nbreaks = 6,
       leg_title = "Quantiles",
       leg_val_rnd = 1)
mf_layout("Y",
          frame=T,
          credits="",
          arrow=F)

# Carte de Ystd
mypal<- brewer.pal(6, "RdYlBu")
mybreaks<-c(-10,-2,-1,0,1,2,10)
mf_map(sel, var="Y_std", type="choro",
       breaks=mybreaks, 
       pal = mypal,
       leg_title = "moy. & ec.type",
       leg_val_rnd = 1
       )
mf_layout("Y_std",
          frame=T,
          credits="",
          arrow=F)

Valeurs de voisinage

sel$Y_lag<-as.numeric(lag.listw(map_nb_w,sel$Y))
sel$Y_std_lag<-as.numeric(lag.listw(map_nb_w,sel$Y_std))
#mapsel$Y_std_lag<-lag.listw(map_nb_w,mapsel$Y_std)
head(st_drop_geometry(sel))
  code                     nom        Y      Y_std    Y_lag  Y_std_lag
1   01                     AIN 36.53036  0.3548758 31.39537 -0.3338932
2   02                   AISNE 50.63966  2.2473937 40.51463  0.8892969
3   03                  ALLIER 36.69576  0.3770614 35.98661  0.2819411
4   04 ALPES-DE-HAUTE-PROVENCE 36.30369  0.3244723 37.09966  0.4312381
5   05            HAUTES-ALPES 31.05594 -0.3794221 33.09124 -0.1064220
6   06         ALPES-MARITIMES 37.73249  0.5161209 39.76664  0.7889666

Indice et diagramme de Moran

spdep::moran.test(sel$Y,map_nb_w, alternative = "two.sided")

    Moran I test under randomisation

data:  sel$Y  
weights: map_nb_w    

Moran I statistic standard deviate = 8.8275, p-value < 2.2e-16
alternative hypothesis: two.sided
sample estimates:
Moran I statistic       Expectation          Variance 
       0.56608333       -0.01075269        0.00427000 
par(mfrow=c(1,1))
moran.plot(x=sel$Y_std,
           listw=map_nb_w, 
           pch=20,
           main ="Diagramme de Moran",
           xlab = "Variable observée (Y)", 
           ylab=  "Variable de voisinage (Ylag)",labels = F)
text(sel$Y_std, 
     sel$Y_std_lag, 
     sel$nom, 
     cex=0.5, 
     col="red",
     pos=1)

Calcul des LISA

locm<-localmoran(sel$Y,map_nb_w,alternative = "two.sided")
sel2<-as.data.frame(locm)
tabres<-cbind(sel,sel2)
head(st_drop_geometry(tabres)) 
  code                     nom        Y      Y_std    Y_lag  Y_std_lag
1   01                     AIN 36.53036  0.3548758 31.39537 -0.3338932
2   02                   AISNE 50.63966  2.2473937 40.51463  0.8892969
3   03                  ALLIER 36.69576  0.3770614 35.98661  0.2819411
4   04 ALPES-DE-HAUTE-PROVENCE 36.30369  0.3244723 37.09966  0.4312381
5   05            HAUTES-ALPES 31.05594 -0.3794221 33.09124 -0.1064220
6   06         ALPES-MARITIMES 37.73249  0.5161209 39.76664  0.7889666
           Ii         E.Ii     Var.Ii       Z.Ii Pr.z....E.Ii..
1 -0.11976472 -0.001368720 0.02025013 -0.8319995     0.40540921
2  2.02009056 -0.054893417 0.76861540  2.3667926     0.01794298
3  0.10745221 -0.001545205 0.02285718  0.7209496     0.47094050
4  0.14142937 -0.001144240 0.02055288  0.9944955     0.31998170
5  0.04081302 -0.001564614 0.03551380  0.2248734     0.82207774
6  0.41158068 -0.002895109 0.13420143  1.1314115     0.25788193

Carte des LISA

par(mfrow=c(1,2))

# Carte de Ystd
mypal<- brewer.pal(6, "RdYlBu")
mybreaks<-c(-10,-2,-1,0,1,2,10)
mf_map(sel, var="Y_std", type="choro",
       breaks=mybreaks, 
       pal = mypal,
       leg_title = "moy. & ec.type",
       leg_val_rnd = 1
       )
mf_layout("Variable standardisée",
          frame=T,
          credits="",
          arrow=F)


mypal<-brewer.pal(5,"BrBG")
mybreaks<-c(-10,-2,-1,1,2,10)
mf_map(tabres, var="Z.Ii", type="choro",
       breaks=mybreaks, 
       pal = mypal,
       leg_title = "Z-score",
       leg_val_rnd = 1
       )
mf_layout("LISA",
          frame=T,
          credits="",
          arrow=F)

Typologie de Moran

# create a new variable identifying the moran plot quadrant for each observation, dismissing the non-significant ones
tabres$typ <- NA

# high-high quadrant
tabres[(tabres$Y_std >= 0 & 
                 tabres$Y_std_lag >= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "high-high"
# low-low quadrant
tabres[(tabres$Y_std <= 0 & 
                 tabres$Y_std_lag <= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "low-low"
# high-low quadrant
tabres[(tabres$Y_std >= 0 & 
                 tabres$Y_std_lag <= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "high-low"
# low-high quadrant
tabres[(tabres$Y_std <= 0 & 
                 tabres$Y_std_lag >= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "low-high"
# non-significant observations
tabres[tabres$Pr.z....E.Ii.. > 0.05, "typ"] <- "not signif."  

tabres$typ <- as.factor(tabres$typ)
tabres$typcol<-NA
tabres$typcol <- case_when(tabres$typ=="high-high" ~ "red",
                        tabres$typ=="low-low" ~ "blue",
                        tabres$typ=="high-low" ~ "pink",
                        tabres$typ=="low-high" ~ "lightblue",
                        tabres$typ=="not signif." ~ "gray90"                      
                        )

mf_map(tabres, type="base",  col=tabres$typcol)
mf_layout("Typologie de Moran", 
          frame=T,
          credits ="")
mf_legend(type="typo",val = c("haut-haut","bas-bas","haut-bas","bas-haut","non sign."),pal = c("red","blue","pink","lightblue","gray90"),title = "Types")

C. ECHELLE DES CIRCONSCRIPTIONS

On passe à l’échelle des circonscriptions législatives qui sont de niveau infradépartemental. On exclue le département de Paris (75) dont les données sont imparfaites.

Préparation des données

### Fonds de carte
map<-st_read("data/elect2024/map_circ.shp", quiet=T) %>%
    filter(code_dpt!=75) %>%
      select(code=ID,  geometry) %>%
      st_transform(2154)

### Données
don<-readRDS("data/elect2024/don_circ.RDS") %>%
          filter(code_dpt!=75) %>%
          select(code=ID, nom=ID, 10:20)


### Jointure
mapdon<-left_join(map[,1],don)

# Table de voisinage
map_nb<-spdep::poly2nb(mapdon,row.names = mapdon$code)

# Table de poids
map_nb_w<-nb2listw(map_nb)
summary(map_nb_w)
Characteristics of weights list object:
Neighbour list object:
Number of regions: 517 
Number of nonzero links: 2776 
Percentage nonzero weights: 1.038576 
Average number of links: 5.369439 
Link number distribution:

  1   2   3   4   5   6   7   8   9  10  11 
  3  18  64  79 116 100  73  46  11   6   1 
3 least connected regions:
17001 56005 59013 with 1 link
1 most connected region:
42004 with 11 links

Weights style: W 
Weights constants summary:
    n     nn  S0       S1       S2
W 517 267289 517 208.7879 2132.706
# Carte de voisinage
par(mfrow=c(1,1))
coo<-st_coordinates(st_centroid(mapdon))
mf_map(mapdon, type="base",col="lightyellow")
mf_layout("Carte des liens de voisinage", frame=T)
plot.nb(map_nb,coords = coo,add = T,col = "red",points = F)

#mf_label(mapdon, var="code", col="blue",cex=0.6,halo = T, bg="white",r = 0.1)

Choix d’une variable Y

sel<-mapdon %>% select(code,nom, Y=varY) %>%
                mutate(Y_std = as.numeric(scale(Y)) ) 
par(mfrow=c(1,2))
# Carte de Y
mf_map(sel, var="Y", type="choro",
       nbreaks = 6,
       lwd=0.1,
       leg_title = "Quantiles",
       leg_val_rnd = 1)
mf_layout("Y",
          frame=T,
          credits="",
          arrow=F)

# Carte de Ystd
mypal<- brewer.pal(6, "RdYlBu")
mybreaks<-c(-10,-2,-1,0,1,2,10)
mf_map(sel, var="Y_std", type="choro",
       breaks=mybreaks, 
       pal = mypal,
       lwd=0.1,
       leg_title = "moy. & ec.type",
       leg_val_rnd = 1
       )
mf_layout("Y_std",
          frame=T,
          credits="",
          arrow=F)

Valeurs de voisinage

sel$Y_lag<-as.numeric(lag.listw(map_nb_w,sel$Y))
sel$Y_std_lag<-as.numeric(lag.listw(map_nb_w,sel$Y_std))
#mapsel$Y_std_lag<-lag.listw(map_nb_w,mapsel$Y_std)
head(st_drop_geometry(sel)) 
   code   nom        Y      Y_std    Y_lag  Y_std_lag
1 33004 33004 30.57356 -0.1781933 28.62401 -0.3751887
2 38001 38001 15.17156 -1.7345108 24.45730 -0.7962199
3 59010 59010 34.87720  0.2566749 25.39481 -0.7014880
4 33007 33007 19.95650 -1.2510096 23.81798 -0.8608204
5 01001 01001 37.79588  0.5515966 37.05787  0.4770237
6 01002 01002 36.73677  0.4445777 33.63850  0.1315078

Indice et diagramme de Moran

spdep::moran.test(sel$Y,map_nb_w, alternative = "two.sided")

    Moran I test under randomisation

data:  sel$Y  
weights: map_nb_w    

Moran I statistic standard deviate = 26.255, p-value < 2.2e-16
alternative hypothesis: two.sided
sample estimates:
Moran I statistic       Expectation          Variance 
     0.7284488328     -0.0019379845      0.0007739103 
par(mfrow=c(1,1))
moran.plot(x=sel$Y_std,
           listw=map_nb_w, 
           pch=20,
           main ="Diagramme de Moran",
           xlab = "Variable observée (Y)", 
           ylab=  "Variable de voisinage (Ylag)",labels = F)
text(sel$Y_std, 
     sel$Y_std_lag, 
     sel$nom, 
     cex=0.5, 
     col="red",
     pos=1)
text(-1,-0.3,"Low-Low",col="grey30")
text(-1,0.3,"Low-High",col="grey30")
text(1,-0.3,"High-Low",col="grey30")
text(1,0.3,"High-High",col="grey30") 

Calcul des LISA

locm<-localmoran(sel$Y,map_nb_w,alternative = "two.sided")
sel2<-as.data.frame(locm)
tabres<-cbind(sel,sel2)
head(st_drop_geometry(tabres))
   code   nom        Y      Y_std    Y_lag  Y_std_lag          Ii          E.Ii
1 33004 33004 30.57356 -0.1781933 28.62401 -0.3751887  0.06698569 -0.0000616558
2 38001 38001 15.17156 -1.7345108 24.45730 -0.7962199  1.38372844 -0.0058417792
3 59010 59010 34.87720  0.2566749 25.39481 -0.7014880 -0.18040333 -0.0001279258
4 33007 33007 19.95650 -1.2510096 23.81798 -0.8608204  1.07898156 -0.0030388720
5 01001 01001 37.79588  0.5515966 37.05787  0.4770237  0.26363456 -0.0005907915
6 01002 01002 36.73677  0.4445777 33.63850  0.1315078  0.05857873 -0.0003837836
       Var.Ii       Z.Ii Pr.z....E.Ii..
1 0.004500391  0.9994395     0.31758186
2 0.996965365  1.3916835     0.16401828
3 0.021957452 -1.2165935     0.22375891
4 0.310831367  1.9407654     0.05228674
5 0.060577565  1.0735409     0.28302851
6 0.028004175  0.3523422     0.72458168

Carte des LISA

par(mfrow=c(1,2))

# Carte de Ystd
mypal<- brewer.pal(6, "RdYlBu")
mybreaks<-c(-10,-2,-1,0,1,2,10)
mf_map(sel, var="Y_std", type="choro",
       breaks=mybreaks, 
       pal = mypal,
       lwd=0.1,
       leg_title = "moy. & ec.type",
       leg_val_rnd = 1
       )
mf_layout("Variable standardisée",
          frame=T,
          credits="",
          arrow=F)


mypal<-brewer.pal(5,"BrBG")
mybreaks<-c(-10,-2,-1,1,2,10)
mf_map(tabres, var="Z.Ii", type="choro",
       breaks=mybreaks, 
       pal = mypal,
       lwd=0.1,
       leg_title = "Z-score",
       leg_val_rnd = 1
       )
mf_layout("LISA",
          frame=T,
          credits="",
          arrow=F)

Typologie de Moran

# create a new variable identifying the moran plot quadrant for each observation, dismissing the non-significant ones
tabres$typ <- NA

# high-high quadrant
tabres[(tabres$Y_std >= 0 & 
                 tabres$Y_std_lag >= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "high-high"
# low-low quadrant
tabres[(tabres$Y_std <= 0 & 
                 tabres$Y_std_lag <= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "low-low"
# high-low quadrant
tabres[(tabres$Y_std >= 0 & 
                 tabres$Y_std_lag <= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "high-low"
# low-high quadrant
tabres[(tabres$Y_std <= 0 & 
                 tabres$Y_std_lag >= 0) & 
                tabres$Pr.z....E.Ii..<0.05, "typ"] <- "low-high"
# non-significant observations
tabres[tabres$Pr.z....E.Ii.. > 0.05, "typ"] <- "not signif."  

tabres$typ <- as.factor(tabres$typ)
tabres$typcol<-NA
tabres$typcol <- case_when(tabres$typ=="high-high" ~ "red",
                        tabres$typ=="low-low" ~ "blue",
                        tabres$typ=="high-low" ~ "pink",
                        tabres$typ=="low-high" ~ "lightblue",
                        tabres$typ=="not signif." ~ "gray90"                      
                        )

mf_map(tabres, type="base",  col=tabres$typcol)
mf_layout("Typologie de Moran", 
          frame=T,
          credits ="")
mf_legend(type="typo",val = c("haut-haut","bas-bas","haut-bas","bas-haut","non sign."),pal = c("red","blue","pink","lightblue","gray90"),title = "Types")